﻿using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;

namespace PowerlawDistribution
{
    /// <summary>
    /// This class will generate random numbers following a power law distribution
    /// </summary>
    public class RandomNumbers
    {
        private Random rnd = new Random(System.Environment.TickCount);
        double imax;
        double v;
        private double q, oneMQ, oneMQ_Inv, hxHxm, hxm, s;

        public RandomNumbers(int imax)
        {
            this.s = 1.0;
            this.imax = imax-1;

            init();
        }

        public RandomNumbers(double s, int imax)
        {
            this.s = s;
            this.imax = imax-1;

            init();
        }

        private void init()
        {
            q = s;
            if (q <= 1.0)
                q = 1.1;
            
            if (v < 1)
                v = 1;            

            oneMQ = 1.0 - q;
            oneMQ_Inv = 1.0 / oneMQ;
            hxm = H(imax + 0.5);
            hxHxm = H(0.5) - Math.Exp(Math.Log(v) * (-q)) - hxm;
            s = 1 - Hinv(H(1.5) - Math.Exp(-q * Math.Log(v + 1.0)));
        }

        private double H(double x)
        {
            return (Math.Exp(oneMQ * Math.Log(v + x)) * oneMQ_Inv);
        }

        private double Hinv(double x)
        {
            return Math.Exp(oneMQ_Inv * Math.Log(oneMQ * x)) - v;
        }

        public long next()
        {
            while (true)
            {
                double ur = hxm + rnd.NextDouble() * (hxHxm);
                double x = Hinv(ur);

                long k = (long)(x + 0.5);
                if (k - x <= s)
                    return k;
                else if (ur >= (H(k + 0.5) - Math.Exp(-Math.Log(k + v) * q)))
                    return k;                
            }
        }
        
       
    }
}
